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Abstract. We present a time-dependent radiative model of the atmosphere of HD209458b 
' and investigate its thermal structure and chemical composition. In a first step, the stellar 

^ ' heating profile and radiative timescales were calculated under planet-averaged insolation 

5— I ' conditions. We find that 99.99% of the incoming stellar flux has been absorbed before reach- 

ing the 7 bar level. Stellar photons cannot therefore penetrate deeply enough to explain the 
large radius of the planet. We derive a radiative time constant which increases with depth 

•rH . 

and reaches about 8 hr at 0.1 bar and 2.3 days at 1 bar. Time-dependent temperature pro- 
files were also calculated, in the limit of a zonal wind that is independent on height (i.e. 
solid-body rotation) and constant absorption coefficients. We predict day-night variations 
of the effective temperature of ~ 600 K, for an equatorial rota tion rate of 1 km s~\ in good 



agreement with the predictions by 



Showman & Guillot 



20021). This rotation rate yields day- 



to-night temperature variations in excess of 600 K above the 0.1-bar level. These variations 
rapidly decrease with depth below the 1-bar level and become negligible below the ~5-bar 
level for rotation rates of at least 0.5 km s -1 . At high altitudes (mbar pressures or less), the 
night temperatures are low enough to allow sodium to condense into NaiS. Synthetic transit 
spectra of the visible Na doublet show a much weaker sodium absorption on the morn- 
ing limb than on the evening limb. The calculated dimming of the sod i um fe ature during 
planetary transites agrees with the value reported by 



Charbonneau et al 



2002) 
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1. Introduction 



The discovery of HD209458b dCharbonneau et al 



(2000); 



Henry et al 



(2000)) allows for the 



first time probing directly the structure of a planet outside our Solar System. Indeed, the fact 
that it transits in front of its star allows both the measurement of its radius and the spectroscopic 
observation of its atmosphere. 

Quantitatively, HD209458 is a GO subgiant, with a mas s Af» = 1.06M o , radius R + = l.\8R Q 



and age = 5.2 Gyr, with unce rtainties of 10% or more (Codv & Sasselov 



orbits the star in 3.524739 days dRobichon & ArenoukoOOl) . at a distance ~ 0.045 AU. Models 
of the star imply that its r adius is R = 92 200-109000 km (about 40% more than Jupiter) for a 



2002). The planet 



mass M = 0.69 + 0.02 Mj (Brown et al.1 2001 



The relati vely large radius of the planet appears difficult to explain using standard evolu 



tion models ( Bodenheimeretal 



nl 120021 



2001 



2003 



Guillot & Showman 



2002; 



BarafFe et al 



2003). 



Guillot & Showman (2002) concluded that, for a realistic model of the outer atmosphere irra- 
diated by the parent star, an extra energy source is needed at deep levels to explain HD209458's 
radius. It could result from the penetration of a small fraction of the stellar flux at pressures of tens 
of bars or from kinetic energy produced from the stella r heating and t ranspo rted downwards to 



the interior region. The same conclusion was reached bv lBaraffe et alJ (120031) from model calcu- 
lations directly coupling the irradiated outer atmosphere and the interior. However, the possibility 
that systematic errors both in the determin ation of the stellar radius and atmospheric properties 



of the planet remains (iBurrows et aikoO 

In order to accurately model the evolution of extrasolar giant planets, and hence gain infor- 
mation on their composition, one has to understand how their atmospheres intercept and reemit 
the stellar irradiation. The problem is especially accute for planets like 5 1 Peg b and HD20 9458b 



(herea fter Pegasides) which are believed to be locked into synchronous rotation llGuillot et al 



1996) due to their proximity to their star. As a consequence the amount of irradiation received on 



the day side and the insolation pattern are unlike what is experienced by any planet in the Solar 
System. 

Several studies have investigated the radiative equilibrium structure of Pegasides, treating the 
atmosphere as a one-dimensional column receiving an average flux from above and a smaller 
intrinsic flux from below. In these models, the stellar heat is either even ly distributed over the 



entire planet (Seager& Sasselov 



1998 



over the day side ( Sudarskv et al 



2003 



2000 



redistributed at all (Barman et al. 2001 



Gouken 



Baraffe et al 



euaue 



2003 



et al 



2000), o r redis tributed only 



Burrows et al 



2003), or even not 



Burrows et 



alfcooJl. 



In fact, in a planetary atmosphere, 



winds carry part of the stellar heat from the day side to the night side and from equator to poles, 
so that the thermal structure depends not only on the insolation pattern but also on the dynamics. 
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Two recent dynamical studies ( Showman & Guillo 



2002; 



Cho et al 



2003) have shown 



that the atmospheric structure probably presents strong latitudinal and longitudinal variations 
in temperature (and hence composition). These two investigations differ in several respects: 



Showman & Guillo used analytical arguments to show that the atmospheres of Pegasides 

should have strong winds (~ 1 kms~') and relatively strong day-night and equator-to-pole tem- 
perature contrasts (~ 500 K near optical depth unity). These estimates are confirmed by prelimi- 
nary 3D simulations done using the EPIC (Explicit Planetary Isentropic Coordinate) model with 
a radiative time constant of 2.3 days. These simulations tend to yield a prograde equatorial jet 
and thus imply that t he atmosph ere could superrotate as a whole (as is the case for Venus). On 



the other hand, 



Cho et al 



( 2003) solved 2D shallow-water equations, assuming a characteristic 
wind speed, related to the mean kinetic energy, of 50 to 1000 ms -1 and a radiative equilibrium 
time of 10 days. Their simulations yield a circulation which is characterized by moving polar 
vortexes around the poles and ~ 3 broad zonal jets. For the highest speed considered (1 km s~'), 
the temperature minimum is about 800 K and the maximum temperature contrast is ~ 1000K. 

The goal of the present study is to apply a time-dependent radiative transfer model to the 
case of HD209458b in order to: (i) determine its mean temperature structure and stellar heating 
profile; (ii) estimate the characteristic radiative heating/cooling timescale as a function of depth; 
(iii) model longitudinal temperature variations assuming a solid rotation; (iv) infer consequences 
for the variations in chemical composition in light of the spectroscopic transit observations. 

In Section |2 we present our time-dependent radiative transfer model. We then apply the 
model to HD209458b. In Section [3] we compare models obtained for averaged insolation con- 
ditions, and calculate the corresponding radiative timescales. In Section |4] we calculate the 
longitude-dependent thermal structure of the atmosphere by assuming a solid body rotation, 
mimicking a uniform zonal wind. The variations in the chemical composition induced by temper- 
ature variations are investigated with the emphasis on the condensation of sodium. A summary 
of the results and a conclusion are presented in Section|5] 



2. The Atmospheric Model 

2.1. Physical problem 

In a one-dimensional model, the evolution of the temperature profile in an atmosphere under 
hydrostatic equilibrium is related to the net flux F by the following energy equation: 

dT^^mgldF] 
dt C p \dpj ' 

where mg is the mean molecular weight and C p is the mean specific heat. The net flux is divided 
into the thermal flux emitted by the atmosphere Fir (upward - downward) and the net stellar flux 
F* (downward - upward), so that F = Fir - F*. Equation^can then be rewritten as: 

*L=(h(p)-c(p)) , (2) 
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where h{p) — —jr-^j*- is the heating rate and c(p) — is the cooling rate. 

Radiative equilibrum corresponds to a steady-state solution of the energy equation and is thus 
obtained by setting the left-hand term of Eqs. 1-2 to zero. In this case, the flux is conservative 
and heating and cooling rates are equal at any level in the atmosphere. In a first step, we calculate 
the radiative equilibrium solution for HD209458b, using a planetary-averaged stellar irradiation. 
In a second step, we investigate the variations in the temperature profile due to the time-varying 
insolation, which requires solving Eq.|2] 



Goukenleuque et al 



( 2000), with, in the present model, 



2.2. Numerical method 

We use the atmosphere code described in 
an atmospheric grid of yV = 96 levels from 3 x 10 3 to 1 x 10~ 6 bar. To calculate F* as a function 
of pressure level p, we solve the radiative equation of transfer with scattering in the two-stream 
approximation and in plane-parallel geometry using a monochromatic line-by-line code. The 
boundary condition is that the incident downward flux F^(0) is given by: 

Fi(0) = «r(^) B Y (T+) , (3) 

where /?* is the star's radius, a the distance of the planet to the star's surface, and B v (r*) is the 
monochromatic Planck function at the star's brightness temperature r*. To solve for steady-state 
radiative equilibrium, we consider disk-averaged insolation conditions, and thus use a — 4. For 
our time-dependent calculations, we use a = n max[cos(/l),0], where A is the zenith angle of the 
star. The stellar flux is calculated from 0.3 to 6 pm (1700-32000 cm -1 ). It is assumed that short- 
ward of 0.3 pm photons are either scattered conservatively or absorbed above the atmospheric 
grid (p < 1 /ibar) and thus do not participate to the energy budget. 

The planetary thermal flux is calculated from 0.7 to 9 pm. At level p\, it can be expressed as: 

N 

FlR = ZZ fli J- kBk ( r j) ' (4) 

* 7=1 

where (r^j is the Planck function at the temperature Tj of pressure level pj and at frequency Vk, 
and flijjt (= Oj^k) is a coupling term relat ed to the transmittance betw een levels p\ and pj, averaged 
over a frequency interval centered at ( Goukenleuque et alJ2000 ). We use an interval width of 
20 cm" . These coefficients are calculated through a line-by-line radiative transfer code with no 
scattering. Below the lower boundary located at level p\, an isothermal layer of infinite optical 
depth is assumed with a temperature T gr . It does not represent the actual temperature profile in 
the interior. In fact, for numerical reasons, what is exactly set in the model is the upward flux at 
the lower boundary. 



2.2.1 . Steady-state case 

Starting from an initial guess temperature profile and the associated abundance gas profiles, 
stellar and thermal fluxes are computed on the atmospheric grid. The heating and cooling rates, 
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h(p\) and c(pi), defined above are then calculated. At each level of the grid, the temperature is 
modified by an amount: 

ATi = eipi) [h(pd - c( Pi )] + rjipO (T^ - IT, + T i+l ) (5) 

The coefficient e(p), analogous to a time step, increases with pressure to account for the increase 
of the radiative time constant with depth and thus to speed up convergence. Its value is about one 
hour in the upper layers and is as long as few days in the deepest layers. 

The second term in the equation introduces some numerical dissipation needed to avoid in- 
stabilities at high pressure levels where the atmospheric layers are optically thick. The coefficient 
rj(p) is set to 1 x 1CT 3 above the 0.3-bar pressure level and reaches 2% at the lower boundary 
(pi). The temperature T gr beneath pressure level p\ is set to T\ + 5T\, where 6T\ is calculated 
to yield F 1R (pi) = Fi nt . The thermal flux is then calculated using the updated temperature profile 
and the procedure is iterated until the condition h(pi) = c{pi) is fulfilled at each level with a 
precision of 5%. 

A convective adjustment is then applied to the solution profile in regions where the radiative 
gradient exceeds the adiabatic value. In these convective regions, the gradient is set to the adia- 
batic value and the energy flux then includes both convective and radiative components ( Eq. 1 
does not hold anymore). 

After convergence, the equilibrium gas profiles are re-calculated using the solution tempera- 
ture profile. The stellar flux deposition and the fly^ coefficients involved in the calculation of the 
thermal flux are then re-computed and a new temperature profile is obtained from the iterative 
process described above. If this new solution profile is close enough to the previous one (with a 
precision of 10 K), we retain it as the radiative equilibrium solution. If not, the whole iterative 
process is continued. 

2.2.2. Time-dependent case 

To investigate the atmospheric response to a time-varying insolation, we set the incident stellar 
flux according to Eq. |3]with a - n max[cos(/l),0] and A = ^f-, T being the rotation period of 
the atmosphere. This insolation pattern represents the variation of the stellar flux during the day, 
averaged over latitude (i.e., along a meridian). It yields the same daily insolation as the steady- 
state planet-average case (< a > t = |). We then solve Eq. fusing a time-marching algorithm 
with a timestep of 300 sec, starting from the radiative-equilibrium solution profile. Integration 
is performed over several rotation periods until the temperature reaches a periodic state at each 
level. 

2.3. Opacities 

We include the following opacity sources: Rayleigh scattering (for the stellar flux component 
only), collision-induced absorption from H2-H2 and H2-He pairs, H bound-free, free-free 
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absorption, molecular rovibrational bands from H2O, CO, CH4, and TiO, and resonance lines 
from Na and K. Spectrosco pic data used for H2-H2, H2-He, H2O, CO, and CH4 are described in 
Goukenleuaue et al.l (12000I) . An improvement over Goukenleuque et al.'s modeling was to add 



absorption from the neutral alkali metals Na and K whose important role was established by 



Burrows et al 



(2000). For the resonance lines at 589 nm (Na) and 770 nm (K), we followed the 



general prescription of 



Burrows et al 



(1200 0) to calculate their collisionally broadened lineshape. 
In the impact region, within Act of line center, a Lorentzian profile is used with a halfwidth 
calculated from the simple impact theory: y = 0.071(T/2000)~ ' 7 cirT 1 ataT 1 for Na and y = 
0.14 (T/2000r a7 cm" 1 atnr 1 for K. 

Beyond the transition detuning Act, we use a lineshape varying as (v - vo) -3 ^ 2 as predicted 
from the statistical theory. We multiply this lineshape by e - h ( v - v oV kT to account for the expo- 
nential cutoff of the profile at large fre quencies (~ 1000-300 cirT 1 at 2000 K). We adopt the 



detuning frequencies Act calculated by 



Burrows et al 



(2000), 30 cnT 1 (7/500 K) for Na and 



20 cm" "'(77500 K) ub for K. H bound-free and fre e-free absorption were modeled as in 



Guill ot et al 



( 1994) using the prescription of 



John 



(1988) for H and 



om lBretJ dl 

have limitations, especially for effective temperatures ~ 2500 K 1 Allard et alJ2000 ). This should 
have a minor impact on the structure of HD209458b because of its relatively low effective tem- 
perature, and appears justified in regard of the other sources of uncertainties of the problem. 



Belli ill 98(1) for H'. 



A relatively simplified chemical equilibrium is calculated using the ATLAS code llKurucz 



19701) including the condensation of Na (as Na2S), K (as K2S), Fe (as Fe), Mg and Si (as 



MgSiOs), and assuming that Al, Ca, Ti and V also condense with simila r pa rtial pressures a s 

t Fe glev & Lodders 1 1994) and lLodderJJl999l) . 



MgSi03. The condensation curves are taken from J 

An important uncertainty in any attempt to calculate the atmospheric structure of Pegasides 
concerns the presence and structure (in term s of composition and grain sizes) of silicate clouds 



in the atmosphere (e.g. 



Fortnev et al 



20031) . As discussed by 



Showman & Guilloj J2002h . such 



a cloud could either occur on the day side or on the night side depending on whether the global 
atmospheric circulation yields a strong vertical advection at the substellar point or a mainly 
superrotating atmosphere. Given this unknown and the very uncertain microphysics underlying 
the cloud structure, we choose not no include any additional opacity due to condensed particules 
in the present model. 



2.4. Input parameters 



1989). The acceleration of 



We assume a solar abundance of the elements (lAnders & Grevessd 
gravity was set to g — 9.7 m s~ 2 correspo nding to a radius of 1.35 Jupiter radius and a mass of 



0.7 time that of Jupiter ( MazdietaL|200C )■ To calculate the incoming stellar flux, we adopt R+ 



Sun given in Pierce & Allen 



1.18 R B , a = 0.047 AU 1 Mazeh e { al.l2000l) . and use the brightness temperature spectrum of the 



N. Iro, B. Bezard and T. Guillot: A Time-Dependent Radiative Model of HD209458b 7 

Our lower boundary condition is fixed by the intrinsic flux of the planet = crT^.. 
Standard evolution models, in which the stellar flux is totally absorbed at low pressure levels 
(< 10 bars), indicate an intrinsic effective temperature T\ M of approximatively 100 K at the age 
of HD209458b. On the other hand, to reproduce the rela tively large observed radius of the planet, 



an int rinsic temperature as high as ~300 K is needed l lGuillot & Showmanl2002 ; 



Baraffe et al 



2003), which requires an extra source of energy in the interior. In our model, we have thus con- 
sidered two possible boundary conditions corresponding to T mt = 100 K (cold case) and 300 K 
(nominal case). 



3. Radiative equilibrium solutions 

3.1. Temperature profiles 

We compute two atmospheric temperature profiles corresponding to the two boundary conditions 
(Tint = 100 K and T mX = 300 K), as shown in Fig. ^ Above the 10-bar level, the boundary 
condition does not make any difference because the incoming stellar flux is dominant with respect 
to the intrinsic heat flux. A convective zone appears below the 1-kbar level for the cold case and 
the 0. 1-kbar level for the nominal case. 

We infer an effective temperature 1 T e g = 1340 and 1350 K for the two boundary conditions 
(cold and nominal respectiv ely), and a planetary Bond albedo of 0.10. While this albedo closely 
agrees with that inferred bv lBaraffe et alJ II2003I) . our T e g is lower because these authors redis- 
tribute the stellar flux only over the day side and thus have a stellar heating twice larger than ours 
Below the 3-mbar level, our profile is ~ 600 K cooler than that produced by 



twice larger tha 


n ours 


Sudarsky et al. 


(2003 



(as shown in their Fig. 26). This difference partly results from the twice larger insolati on used by 

(2003[KHowever 



Baraffe et al 



these authors who chose the same redistribution of stellar heat as 
it can only produce 15-20% larger temperatures as illustrated by iSudarskv et alJ ( 120031) ' s com- 
parison of profiles with different incident flux weighting (shown in their Fig. 28). The expected 
temperature increase below the 3-mbar level is thus limited to 150-300 K. On the other hand, 
we note that these authors used an intrinsic temperature of the planet of 500 K, hotter than ours, 
which can contribute to the difference with our model at least below the ~0.5-bar level. Above 
the 0. 3-mbar level, the difference is only ~150 K and can be accounted for by the factor of 
two in the stellar flux. Another difference be tween the two models is the incorporation of sili- 
cate and iron clouds in 



Sudarskv et al 



l|2003j)'s model whereas ours is cloud-free. These clouds, 
with bases at 5-10 mbar in their model, result in a cooler atmosphere below the ~30-mbar level 
and a hotter atmosphere above, and therefore are n ot the source of the d iscrepancy in the lower 
atmosphere. The (cloud-free) profile calculated bv iBaraffe et alJ ( 120031) is quasi-isothermal at 
-1800 K between 10~ 4 and 100 bar whereas ours increases from 1000 to 1800 K for the same 



7int=100 K. Lowering the insolation by half in 



Baraf fe et al 



( 2003 )'s model shou ld bring it down 



to ~1500 K but cannot produce the gradient we have or the even steeper one of 



Sudarskv et al 



1 defined from the total flux emitted by the planet by: Fir = <tT* s 



8 
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1000 2000 3000 4000 
Temperoture (K) 

Fig. 1. Temperature profiles of the static disk-averaged solutions for the two boundary conditions. 
There is no difference above the 10-bar level as the intrinsic heat flux becomes negligible with 
respect to the net stellar flux. 



2003 



Be low the 0.3-bar regi on, our profile is intermediate between those of 



2003) and 



Baraffe et al 



Sudarskv et al 



(2003) after correcting to first order for their twice larger stellar heat- 



ing. 



The preferred boundary condition used by 



Guill ot & Showman! 00021) for HD209458b's ra 



dius is: T(p = 3bar) = ri so i a ted(7 , e ff,g) - 1000 K, where ri so i ate d(r e ff, g) is the temperature of an 
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isolate d planet of same effective temperature T e ff and gravity g as calculated by 



Burrows et al 



( 1997). This is relatively close to our solution profile, which implies a 3-bar temperature that is 



1200 K less than in the isolated case. 



3.2. Penetration of the stellar flux and spectra 



Figure |2 shows the net stellar flux as a function of pressure level for the "cold" and "nominal" 
profiles. We find that 90% of it is absorbed at the ~0.9-bar level, 99% around the 2-bar level, 
and 99.99% at the ~7-bar (resp. 5-bar) level in the cold (resp. nominal) case. At this level and 
deeper, the stellar flux becomes smaller than the intrinsic flux (calculated in the cold case) and 
cannot affect the temperature profile in a significant way. At deep levels, the stellar flux decreases 
more rapidly with depth in the nominal case than in t he col d case due to the larger abundances 



of H 2 and H ions and of TiO. 



Guillot & Showman 



(2002) estimated that a penetration of 1% 



of the stellar flux to p ~ 100 bar in the "cold" case would allow the radius of HD209458b 
to be explained without any other energy dissipation. Clearly our calculations indicate that this 
is not the case and that the large atmospheric opacity beneath pressure levels of a few bars 
prevents any significant fraction of the stellar heat from reaching directly the ~100 bar region. 
Therefore, provided that the error bars on the planet's radius are not underestimated, we confirm 



the need for an addit i onal h eat s ource at deep levels as advocated by 



Guillot & Showm an (2002P and 



Baraf fe et al 



(2003) 



Bodenheimer et al 



(2001) 



The reflected and thermal spectra of the planet are shown in Fig. [5] The thermal spectrum 
is dominated by the water vapor bands, but CO absorption bands around 2100, 4300, and 6400 
cirT 1 are visible. It exhibits spectral windows centered at 2600 crrT 1 (3.8 pm), 4500 cirT 1 
(2.2 pm), 6000 cnr 1 (1.7 pm), 7800 cnr 1 (1.28 pm), and 9400 cnr 1 (1.07 pm). The most 
intense one at 3.8 pm is limited on the short-frequency side by the (1-0) CO band and on the 
high-frequency side by the V3 H2O band. The stellar reflected spectrum peaks around 28 500 
crrT 1 (0.35 pm). At wavenumbers less than ~18 000 cm -1 , it is almost completely absorbed due 
to weaker Rayleigh scattering and strong atmospheric absorption. In particular, in the region of 
the alkali lines (13 000-17 000 cm -1 ), the flux emitted by the planet is at minimum while that of 
the star is large. This emphasizes the i mportance of a l kali li ne absorption in the energy balance 



of the Pegasides, as first recognized by 



Burrows et al 



(2000). 



3.3. Radiative Timescales 

To characterize the radiative time constant r lal j at a given pressure level po, we applied a gaussian 
perturbation to the radiative equilibrium temperature profile and let it relax to its equilibrium 
state using Eq. |3 The perturbation has the form AT(p) = A7q 2 I 21 "™] , i.e. a full width at 
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Fig. 2. Net stellar flux as a function of pressure level for the Cold and Nominal profiles. 



half maximum of one scale height. The radiative time constant is calculated from the Newtonian 
cooling equation: 



Ar _ d(AT) 

Trad 



dt 



(6) 



where Ar is the temperature deviation from the equilibrium profile at level po. We verified that 
the Newtonian cooling approximation is justified and yields the same T ra d as far as A7o is small 
(typically less than 5% of the equilibrium temperature). 
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— Thermol emission spectrum 

— Incident stellar spectrum 
■■■■ Reflected stellar spectrum 



5.0x1 3 1.0x10 4 1.5X10 4 2.0x1 4 2.5x1 4 3.0x1 4 
Wavenumber (cm"') 



Fig. 3. Reflected and thermal emission spectra of the static disk-averaged solution for the nominal 
boundary condition, in logarithmic scale. The incident stellar spectrum is shown for comparison. 



The result is shown in Fig. 0] As expected, the radiative time constant increases monot- 
ically with pressure. Arou nd 1 bar, it is ~2.3 days, about twice the value estimated by 



Showman & Guilton ( 120021) . On the other hand, the long radiative timescale (10 days) assumed 



by 



Cho et al 



( 2003) implies that they significantly underestimated the cooling in their circulation 



model. 

Above the 1-bar pressure level, the radiative timescale is relatively short - less than the 
rotation period - and we thus expect a rapid response from the planet's atmosphere to a possible 
atmospheric circulation and large horizontal variations of temperature. In particular, at 0.1 bar, 
T ra d is 3xl0 4 sec (8 hr). This region, which corresponds to the optical l imb of the planet at 0.6 um 
due to Rayleigh scattering, is that probed by the transit observations ( Charbonneau et al.l2002 ). 
The small radiative time constant suggests that the day-night thermal contrast is large and that the 
atmospheric morning limb may be significantly colder than the evening one for any reasonable 
zonal wind speed - as measured in the synchronously-rotating frame -. This likely asymmetry, 
in temperature and hence potentially in the chemical composition, should be kept in mind when 
analyzing spectroscopic transits. Below the 1-bar level, the opacity is large and the timescale 
increases rapidly, varying almost as p 2 . We limited this calculation to 10 bar because below this 
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level our profile becomes convective in the nominal case (Fig.^ and the radiative timescale then 
cannot be calculated with the same method. 



4. Time-dependent solutions for an atmosphere in solid rotation 

4.1. Longitude-dependent temperature profiles 

We then introduce a solid body rotation by moving the atmosphere with a constant angular ve- 
locity with respect to the synchronously-rotating frame. This procedure mimicks a zonal atmo- 
spheric circulation and allows us to investigate the temperature as a function of longitude. The 
effect of rotation is considered only through the modulation of the incoming stellar flux. The 
insolation is maximum at the substellar longitude (noon) and constantly decreases until the at- 
mosphere no longer receives the stellar flux on the night side. This periodic insolation is shown 
in Fig. 13 A first limitation of our model is that heating and cooling rates are valid for layers 
which are static with respect to each other. The rotation of the atmosphere is thus approximated 
by a solid body rotation, i.e. with a zonal wind constant with height. A second limitation is that 
these heating and cooling rates correspond to a given chemical composition (associated with the 
nominal temperature profile in Fig.[0 that does not change during the integration process. 

We compute the temperature profiles for rotation periods corresponding to equatorial wind 
velocities of 0.5, 1 and 2 km sec -1 , in the range predicted bv IShowman & Gui llot (2002J). The 
effective temperature of the planet is shown in Fig. [7] and the temperature at selected levels is 
displayed in Fig.[3]as a function of longitude. Figure[6]shows an equatorial cut of the atmosphere, 
with the color scale for temperature indicated on the right. The temperature maximum shifts with 
respect to the maximum of insolation as depth increases. This phase lag results from the increase 
in the radiative timescale with depth. It is particularly noticeable at p > 1 bar, where T ra( j becomes 
comparable or exceeds the rotation period. We found day/night temperature variations exceeding 
400, 600 and 800 K above the ~0.1-bar level for equatorial win ds of 2, 1 and 0.5 km sec 



Showman & Guillot (2002) 



respectively. These values are of the same order as predicted by 
However, at their reference level of 1 bar and for a wind of 1 km sec 1 , the contrast we get is 
at most 100 K, less than the ~500 K they predict. Part of the discrepancy comes from the half 
lower radiative timescale used by these authors. In fact, the temperature contrast estimated by 
those authors pertains to the level where T ~ T e g which takes place at 0.15 bar in our model, 
instead of 1 bar. Figure 5 shows that we do get contrasts of 400-700 K in T e g for winds in the 
range 2-0.5 km sec -1 . 

Below the ~5-bar level, the temperature field is essentially uniform with longitude (and time) 
due to the long radiative timescale. It is then equal to that given by the radiative equilibrium 
solution under planet-averaged insolation conditions. This points to the need of considering this 
temperature as a boundary condition for evolutionary models rather than t hat cor respondi ng to 



the illuminated hem isphere with no day-to-night redistribution of heat (e.g. . Burrows et al 
BaraffeetalLomi. 



2003 
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An important effect is that the temperature profile becomes cold enough on the night side 
to allow sodium to condense, as shown in Figs.[3H6] The condensation level extends relatively 
deep in the atmosphere on the night side and on the morning limb, down to 0.1, 0.3 and 0.5 bar 
for equatorial wind speeds of 2, 1 and 0.5 km sec -1 respectively. This is close or even below the 
optical limb around the Na resonance doublet at 589.3 nm. Hence, during a transit, we expect 
stellar radiation at this wavelength to be less absorbed on the morning limb than on the evening 
limb. 



4.2. Sodium condensation and transit spectra 



Charbonneau et al 



(2002) conducted HST spectroscopic observations of HD209458 centered 
on the sodium resonance doublet at 589.3 nm. Using four planetary transits, they found a visi- 
ble deeper dimming in a bandpass centered on this feature than in adjacent bands. They inter- 
preted it as absorption from sodium in the planet's atmosphere. However, the sodium absorp- 
tion seen in the data i s lower than p redicted by existing models by approximately a factor three 



(Hubbard etalJl2001 



I2OOI ; Brownb 



Brown 2001). These authors invoke several possibilities to explain such a 



depletion: 



- condensation of sodium into mostly Na2S. However according to their model a very large 
fraction (99%) of the sodium should condense out to explain by itself such a small absorption. 

- Ionization of sodium by the large stellar flux incident on the planetary atmosphere. This 
is a secondary effe ct which cannot explain alone the weak absorption, as confirmed by 



Fortnev et al 



(2003). 



- Very high (above the 0.4-mbar level) particulate opacity. This possibility is difficult to assess 
as modeling of photochemical hazes still needs to be performed. 

- A bulk depletion of the sodium abundance in HD209 458. It appears howeve r that parent stars 



of close-in planets generally have a high metallicity llGonzalezll9S 



2000). 



An alternate explanation is provided by 



Barman et al 



(2002). According to their work, Na 



could be far from being in local thermodynamic equilibrium in the limb region. This effect would 
then reduce the depth of the Na absorption seen during a transit. The amplitude of the effect is 
however difficult to assess due to the lack of well determined collisional deactivation rates from 
H2 and He. 

It should be noted that the models used to an alyze the observations were static and horizon- 



tally uniform whereas 



Charbonneau et al 



(2002)'s data pertain to the planetary limb. Horizontal 
temperature variations due to insolation and dynamics can drive compositional variations that 
should be taken into account to calculate the sodium abundance over the limb. More specifically, 
if sodium condenses on the night side, half of the limb is depleted in sodium above some level in 
the atmosphere. This effect can potentially reduce by half the sodium absorption during a transit 
compared to the prediction from static models. 
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We present here calculations of the spectral variation of the planetary radius in the region of 
the Na lines. The radius r{v) is defined here as the limb corresponding to a tangential extinction 
optical depth of unity. To determine it, we calculate, as a function of frequency, the transmission 
of the atmosphere for a series of lines of sight at radii r;(v) and then interpolate in this grid 
to find r(v). This spectrum r(v) is finally convolved to a resolution of 1000. Two models are 
considered: one, corresponding to the evening limb, incorporates the temperature calculated at 
a phase of n/2 (see Fig. [5} and the associated gas profiles; the second one corresponds to the 
morning limb at a phase of 3n/2, Results are shown in Fig. [S] for an equatorial wind speed of 
1 km sec -1 . As we had anticipated, the Na absorption is much less pronounced in the morning 
spectrum than in the evening one due to sodium condensation. Integrating the spectra over narrow 
ban ds at (588.7-589.9 nm) an d around (581.8-588.7 nm and 589.9-596.8 nm) the Na feature as 
did ICharbonneau et alJ J2002B . we find that the absorption depth is three times weaker on the 
morning limb (6R m = 4.4 x 10 3 7?) compared with the evening limb (6R e = 1.4 x 10 2 R). During 
the transit, the stellar spectrum is absorbed by the whole limb (morning and evening) and our 
model thus predicts a dimming in the sodium band almost half lower than without consensation 
of Na. We calculate it as (6R m + 6R e )R/Rl) = 2.8xl0~ 4 within the uncertainty range of the 
observed value (2.32+0.57xl0~ 4 ). 



5. Conclusions 



We have developped a time-dependent one-dimensional radiative model and applied it to the 
case of HD209458b. Because it is one-dimensional, the model is necessarily simplified in the 
sense that advection is included only horizontally in the limit of a solid-body rotation of the 
atmosphere and vertically with convective adjustements of the temperature profile where appli- 
cable. Furthermore, important effects such as the possible presence of clouds and time-dependent 
variations of the absorption coefficients due to variations of the chemical composition are not in- 
cluded. These effects can be considered in a useful way only through a detailed global circulation 
model, which is beyond the scope of this paper. Conversely, the results of our model can be used 
as a basis for future dynamical models of the atmosphere of HD209458b and of Pegasides. They 
also impact evolution calculations and our perception of the vertical and longitudinal structure of 
strongly irradiated planets, as detailed hereafter. 

First, we point out that when considering constant p lanetary-averaged con ditio ns, our radia 



tive eq uilibrium profile is colder than those calculated by 



Sudarskv et al 



1 20031) and 



Baraffe et al 



(2003) mostly because those authors used a twic e larger insolation. Aft er correcting to first order 



for this difference, our profile agrees wi t h that of Sudarsky^^l 



gion, is still colder than 



Sudarskv et al 



s and lBaraffeetall lE 



(2003) above the ~0.3-mbar re- 



bar (with a maximum discrepancy of 300-400 K), and is intermediate between these two other 
models below this region. These remaining differences probably lie in the opacity sources or the 
radiative transfer treatment but they should not alter the results that follow. 
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In our model, 99.99% of the stellar flux is absorbed at the 5-bar level (99 % at the 2-bar 
level). This rapid absorption prevents a significant fraction (~1 %) of the stellar energy from 
directly reaching pressure levels of several tens of bars. This mechanism c annot thus supply the 



energy source required a t deep levels to explain the radius of HD209458 ( Guil lot & Showman 



2002; 



Baraffeetal. 2003) 



The radiative timescale s that we derive are generally twice longer than tho se used by 



Showman & Gui llot ( 2002) and 5 times shorter than assumed by 



Cho et al 



(2003). Above the 



1-bar level, the radiative time constant is shorter than the rotation period of the planet. These 
relatively short values imply that the atmosphere reacts relatively quickly to perturbations from 
atmospheric dynamics. 

In order to qualitatively test the effect of dynamics on the atmospheric structure, we calculated 
models in which the stellar flux was modulated with a period between 3.5 and 14 days, mimicking 
the effect of a 0.5- to 2-km s _1 equatorial zonal jet on an otherwise synchronously rotating 
planet. 

Depending on the imposed wind, we found longitudinal temperature variations to be between 
400 and 600 K at 0.1 bar, 30-200 K at 1 bar and less than 5 K at 10 bar. This is generally 



con sistent with the re sults obtained by 



that 



tow et. al. 



Showman & Guilloa ( 120021) . On the other hand, the fact 



(2003) obtain temperatures on the night side that can be hotter than on the day 
side is difficult to explain in light of our results. This cannot be ruled out, but would presumably 
imply a combination of strong meridional and vertical circulation. 

The fact that the temperature rapidly becomes uniform with increasing depth implies that the 
mixing most probably takes place in a relatively shallow layer of the atmosphere. In our simula- 
tions, the temperature reached at deep levels is consistent with that of an atmosphere receiving a 
stellar flux averaged over the whole planet on both the day side and the night side (i.e. 1/4 of the 
stellar constant at the planet). This is very important for evolution models, as temperature inho- 
mogeneities at deep levels would t end to fasten the cooling co mpared to homogeneous models 



models (Baraffe et alJ2003 



Burrows et al 



J2003ll 



rtanl^OC 



are calculated from atmospheric boundary con- 



ditions that are relevant to the day side only. These models probably overestimate the temperature 
of the deep atmosphere and therefore the planet's cooling time. 

Finally, the large longitudinal temperature contrast implies that species such as sodium will 
condense on the night side. Even with no settling, the morning limb (which is coldest) can be 
strongly depleted in the condensible species. From calculations of transit spectra representative 
of the morning and evening limbs, we found that the former shows a 3 times weaker sodium 
absorption than the latter. The Na dimming we c alculated through the entir e limb is then in 
agreement with the sodium absorption observed by 
transits. 



Charbonneau et al 



( 2002) during planetary 
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Fig. 4. Radiative timescale as a fu nction of pressure level. Th e symbols show the values of the 
radiative timesca le at 1 bar from Showman & Guillotl l|2002) (1.1 day; diamond) and that of 
Cho etalJ(l2003l) (10 days; square). 
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Fig. 5. Temperature at selected levels as a function of longitude for an equatorial wind velocity 
of (a) 0.5 km sec -1 (b) 1 km sec -1 and (c) 2 km sec -1 . The incoming insolation is shown in 
arbitrary units. Note the shift of the maximum of temperature with respect to the maximum of 
insolation (at the substellar point) as depth (and so the radiative timescale) increases. Below the 
10-bar level, the temperature field is essentially uniform. 
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Fig. 6. Equatorial cut of the atmosphere between the lxl0~ 6 and 10-bar levels for an equatorial 
wind velocity of (a) 0.5 km sec -1 , (b) 1 km sec -1 , and (c) 2 km sec -1 . The level where conden- 
sation of sodium occurs (black line) goes deeper as the night wears on (anti-clockwise) and is 
deepest at the morning limb . Below 10 bar, the temperature field (not shown here) is uniform 
and depends only on the bottom boundary condition. 
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Fig. 7. Effective temperatures as a function of longitude for the three wind velocities. 
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Fig. 8. Spectra of the Na lines region for the morning limb (lower) and the evening limb (upper) 
for a wind velocity of 1 km/s. 



